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Liquid  Rocket  Engine  Combustion  Instability 


FA9550-08-C-0014 


Flow  instabilities  in  roeket  engines  are  triggered  by  the  coupling  between  turbulence, 
flame  fronts  and  aeoustie  waves.  The  ability  to  resolve  seales  of  turbulenee  which  affeet 
the  instability,  tracking  acoustic  wave  motions,  capturing  flame  fronts  and  their 
interaction  with  turbulent  eddies  arc  required  to  numerically  predict  the  onset  of  the 
instability  and  its  limit-cycle  amplitude  and  frequency.  This  report  summarizes  the 
developments  at  Metacomp  Technologies  geared  towards  rocket  engine  flow  simulation. 
The  focus  of  the  eurrent  effort  and  how  it  fits  in  the  overall  goal  of  developing  an 
accurate  unsteady  flow  simulation  capability  is  diseussed  below.  Summary  of  the  effort  at 
Mctacomp  is  followed  by  the  research  at  Penn  State  and  Georgia  Tech. 

Efforts  and  Accomplishments  at  Metacomp  Technologies,  Inc. 

Resolution  of  acoustic  waves  in  low-speed  flows 

Modem  CFD  methods  are  built  on  upwinding  the  inertia  terms  in  the  equations. 
Upwinding  provides  a  systematic  way  of  adding  just  enough  numerical  diffusion  to 
prevent  oscillations  in  the  flow  variables.  Upwinding  requires  knowledge  of  the  local 
waves  speeds.  When  all  the  waves  in  the  system  are  accounted  for,  numerical  diffusion 
couples  all  the  variables.  This  is  also  referred  to  as  matrix  diffusion.  The  main  advantage 
of  this  strategy  is  the  lack  of  tuning  parameters.  However  two  undesirable  effects  show 
up  in  low  speed  flow  simulation.  Firstly,  the  acoustic  speed  which  is  much  large  than  the 
flow  speed  dictates  the  size  of  the  time  step  that  can  be  taken  within  an  explicit  time 
update.  Sinee  the  flow  velocity  is  much  smaller,  a  large  number  of  mueh  smaller  time 
steps  w'ill  have  to  be  computed  before  any  significant  changes  in  the  flow  field  can  be 
deteeted.  It  is  possible  to  overcome  this  with  an  implicit  time  update  procedure.  The 
second  more  detrimental  effeet  arises  from  matrix  diffusion  proportional  to  the  aeoustie 
speed  which  dwarfs  the  flow  speed  on  some  components  of  the  flow.  The  overwhelming 
diffusion  contaminates  the  momentum  equations.  This  translates  to  a  fourth  difference 
smoothing  proportional  to  the  acoustic  speed  on  the  velocity  field. 

Preconditioning  methods  were  developed  to  alleviate  these  two  difficulties.  Artificially 
lowering  the  speed  of  sound  in  the  medium  allows  for  a  larger  time  step  while 
simultaneously  reducing  the  numerical  diffusion  that  affects  low  speed  flow  simulation 
with  the  conservation  form  of  the  equations.  While  preconditioning  works  well  for 
steady  flow  simulation,  a  dual  time  stepping  scheme  is  necessary  for  transient  flow 
simulation  where  tracking  acoustic  waves  in  time  is  eritieal  to  estimating  flow  dynamics 
as  in  a  roeket  engine.  In  the  dual  time  stepping  scheme,  the  conservation  equations 


Q  and  P  refer  to  the  conservation  and  primitive  variables,  respectively,  r  js  the 
preconditioning  matrix.  Convergence  in  pseudo  time  r  satisfies  the  unsteady  form  of  the 
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conservation  equations.  This  approach  is  consistent  with  acoustic  propagation.  However 
in  practice  the  acoustic  waves  are  heavily  damped  even  with  the  dual  time  stepping 
strategy.  A  new  limiter  was  developed  to  control  pressure  damping  which  is  based  on 
local  mesh  sizing.  If  the  local  mesh  size  is  too  large  for  a  given  time  step  size  the 
preconditioning  parameter  is  limited  to  prevent  excessive  smoothing  of  pressure  waves. 

The  following  figure  shows  the  geometry  of  the  test  ease  used  to  assess  the  ability  of  the 
methodology  to  predict  acoustic  modes  in  a  pipe  open  at  the  left  end.  An  oscillating 
boundary  at  the  open  end  generates  standing  waves  within  the  pipe. 


Acoustic  time  for  one  pass  =  0.003  sec. 


Closed 

end 


The  following  figure  shows  the  second  harmonic  of  the  quarter  wave  mode  at  different 
time  instances.  Excellent  agreement  is  obtained  with  analytical  solutions. 


XT  Plar 
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Acoustic  Modes  in  AFRL  Combustion  Stability  Test  Rig 

Simulations  were  also  conducted  to  study  the  acoustic  modes  in  the  AFRL  Combustion 
Stability  test  rig.  The  acoustic  modes  in  the  waveguides  were  estimated.  The  following 
figure  shows  the  pressure  wave  at  130  Hz. 

AFRL  Actuator  -  Mode  1 


» 

The  mode  shape  at  2970  Hz  is  shown  below. 


In  the  two  cases  shown  above,  the  optical  windows  remained  in  place.  The  acoustic 
pathway  is  modeled  as  a  continuous  channel.  The  optical  windows  mesh  seamlessly  with 
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the  walls  of  the  waveguides  at  the  center  of  the  test  seetion.  The  two  drivers  at  either  end 
can  operate  at  syne  or  at  a  specified  phase  difference. 

The  effect  of  close  proximity  of  the  optical  windows  to  the  combustion  zone  was 
addressed  before  hot  firing  the  rig.  Avoiding  the  optical  sections  of  the  wave  guides 
would  be  advantageous  from  the  point  of  view  of  thermal  management.  Removing  the 
optical  seetion  all  together  would  be  desirable  but  its  impaet  on  flow  evolution  was  of 
concern.  The  changes  in  mode  shapes  without  the  optical  sections  were  investigated  by 
modeling  the  entire  volume  of  the  test  seetion. 


The  above  figure  shows  pressure  waves  at  a  frequency  of  2500  Hz.  Modes  from  130  Hz 
to  3000  Hz  were  obtained  to  study  the  chambers  response.  Quieseent  condition  at  room 
temperature  was  assumed  in  evaluating  these  modes. 

Hybrid  RANS/LES  (LNS)  Methodology  for  Unsteady  Flows 

The  LNS  methodology  developed  at  Metaeomp  Technologies  provides  a  systematic  way 
of  LES.  The  LNS  approach  minimizes  sub-grid  dissipation  when  the  computational 
mesh  gets  coarser.  The  LNS  model  takes  mean  flow  quantities  such  turbulence  kinetic 
energy  and  its  dissipation  and  converts  them  into  a  fluctuating  field  of  appropriate 
frequency  when  the  computation  grid  is  sufficient  to  resolve  these  scales.  It  does  not 
need  additional  quantities  whieh  are  difficult  to  measure  at  the  boundaries  to  start  the 
calculation.  When  the  grid  is  sufficiently  fine  all  fluctuation  generated  by  the  velocity 
field  are  resolved.  They  are  converted  to  mean  flow  turbulenee  quantities  only  when  the 
mesh  is  coarse. 
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Temporal  integration 

Preconditioning  is  eoupled  to  a  dual  time  stepping  seheme  to  preserve  conservation  in 
real  time.  The  preeonditioner  operates  on  the  inner  loop  driving  the  solution  to  the  next 
time  level.  Sufficient  number  of  inner  iterations  (<  10)  is  used  to  maintain  second-order 
time  aeeuraey. 

Application  to  a  gas-gas  combustor 

Under  a  separate  effort  funded  by  AFRL,  the  newly  developed  algorithm  diseussed 
above,  was  applied  to  a  2D  eross-seetion  of  a  CH4-O2  gas-gas  atmospherie  pressure  test 
rig1  to  demonstrate  the  ability  to  capture  eombustion  drive  pressure  oscillations.  The 
results  obtained  are  related  to  the  methods  described  above  and  are  presented  here  for 
completeness.  They  demonstrate  the  application  of  the  methodology  to  simulate 
unsteady  reacting  flows. 

The  inlet  geometry  Was  truncated  to  retain  geometric  simplicity.  A  simple  one  step 
chemical  kinetic  mechanism  was  used  in  the  simulation  since  the  goal  was  to  estimate 
whether  the  reecnt  developments  to  CFD++  eould  overcome  exeessive  pressure  damping 
of  preconditioned  systems. 

The  following  figures  showing  the  flow  field  in  the  CH4-O;  burner.  The  figure  on  the  left 
is  with  eombustion  on  and  shows  temperature  eontours  at  an  instant  of  time.  These 
contours  are  raised  to  display  pressure  waves.  When  the  reaction  was  turned  off,  the 
pressure  oscillations  get  quieter  as  shown  on  the  figure  to  the  right.  With  the  eombustion 
off  nearly  eireular  pressure  waves  emanate  from  the  injector  corresponding  to  the 
frequency  at  whieh  vortices  are  shed  off.  CH4  eontours  are  used  in  this  figure  since 
eombustion  is  off  and  temperature  remains  uniform  unlike  the  reacting  flow  shown  on  the 
left. 


Cavitt,  R.C.  ct.  al.,  “On  The  Laboratory  Scale  Survey  of  Pentad  Injector  Stability 
Characteristics”,  A1AA  2007-5587 
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The  output  from  pressure  probe  inside  the  combustor  is  shown  below. 


Pressure  probe  time  history 

Com  bus!  ion  turned  off  at  0.04 1  see 

i - « - 1 - r- - 1  — r — -  •  * — r 

—  low  pas* 


■iL_j _ I _ t _ I _ I _ J _ L- .  -  .  1 _ i _ L_ 

0.02  0.04  0,06  QOS  0.1 

Tunc  \s) 


Combustion  was  turned  off  at  about  0.041  seconds  and  the  green  curve  in  the  above 
figure  shows  the  high  frequency  content  of  the  pressure  tap  decay  in  magnitude  without 
the  heat  release  driving  the  oscillations.  The  low  frequency  oscillations  (shown  in  red) 
mark  bulk  sloshing  within  the  chamber  at  about  70  Hz  after  the  reaction  was  turned  off. 
Vortex  shedding  frequency  is  around  20,  000  Hz  with  significant  amplitude  of  oscillation 
when  reaction  was  on.  The  reaction  also  drives  lower  frequencies  as  shown  by  the  red 
curve  prior  to  0.041  seconds. 
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Heat  Load  on  Wave  Guides 

The  combustor  stability  test  rig  at  AFRL  uses  acoustic  waves  to  perturb  the  fuel  jet  and 
study  its  response  to  acoustic  excitation.  Wave  guides  are  used  to  channel  the  acoustic 
energy  on  to  the  shear  layers  formed  betw  een  the  fuel  and  oxidizer  jets.  Flow  diagnostics 
is  facilitated  by  using  optical  windows  for  the  waveguide  at  the  center  of  the  test  section. 
A  schematic  of  the  facility  is  shown  below. 


CFD++  was  used  to  study  to  simulate  the  flow  through  the  test  section.  A  top  view  of  the 
flow  ehamber  is  shown  in  the  figure  below.  Only  one  quarter  of  the  domain  shown  in  red 
was  modeled  for  computational  efficiency  since  the  steady  temperature  distribution  on 
the  wall  marked  “insulated”  in  the  figure  would  be  a  reasonable  estimate  for  the  heat  load 
on  the  optical  windows. 


CFD++  was  used  to  simulate  steady  flow  in  the  AFRL  test  rig.  The  test  facility  contains 
three  LOX  jets  surrounded  by  annular  gaseous  H2  jets.  The  following  figure  shows  a 
eross-seetion  of  the  combustor.  The  computational  domain  is  shown  in  red. 


Curtain  layer 
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The  mass  flow  rate  was  specified  for  the  LOX  and  gaseous  H2  jets.  This  was  run  at 
stoichiometric  conditions.  Coolant  flow  and  curtain  layer  are  both  employed  as  the  inlet 
boundary  conditions  with  N2  streaming  at  4.5  m/see  at  room  temperature.  In  the 
computational  model  only  three  quarters  of  the  fuel  and  oxidizer  jets  were  simulated. 
LOX  mass  flow  rate  was  1.44  gm/sec  at  100K.  Gaseous  Hydrogen  streams  in  from  the 
annular  jets  at  0.18  gm/see.  Temperature  of  H2  at  inlet  was  200K.  The  LOX  jet  was 
modeled  using  EDP  (Eulerian  Dispersed  Phase)  mode  in  CFD++. 


Information  Sot  Type:  jEquallon  sat  definition 

This  Is  Information  sot  number  1 

The  title  tag  for  this  information  set  Is  eqnsetdefine 

Base  Type. 

Preconditioned  Compressible  RG  NS/Euler  J 

Equation  of  State: 

ideal  Gas  —  [ 

Advection  Terms: 

d 

Viscous  Terms: 

yes  -| 

Turbulence  Simulation: 

♦  RANS  v  LES 

Turbulence  Model: 

2-equation  Realizable  k-eps  —  |  j  More  Models 

Number  of  Color  Tracers: 

1^" 

Number  of  Species  Eqns.  (n-1 ): 

P - 

Fu’ertan  Dispersed  Rase  (EDP): 

♦  yes  v  no 

Number  of  Dispersed  Species: 

li- 

Energy  &  Melting  fraction  for  EDP: 

v  No  ♦  Energy  Only  v  Energy  ♦  Melting  Frac. 

Turn  on  Dtffusfvlty  for  EDP: 

v  yes  ♦  no 

Dispersed  Phase  Model  Type: 

General  EDP  mode  —  j 

Activate  Level  Sets: 

^  ♦ 

Activate  P!  Radiation  Model: 

v  yes  ♦  no 

Conjugate  Heat  Transfer: 

v  yes  ♦  no 

Turbulence  model  recommendations,  for  transonic  Hows  the  Rl  model  or  the 
cubic  k~£  closure  should  be  used.  For  hypersonic  Rows  the  Rt  model  should  be 
invoked.  For  externa!  aero  (where  the  body  Is  far  from  the  in  how  boundary) 
and  for  Internal  flows  the  k-e  R  model  should  be  used.  For  al  other  flows,  incl. 
free  shear  and  free  convection  flows,  the  k  e  or  cubic  k-e  closure  should  be 
Invoked.  FOR  UNSTEADY  FLOWS  USE  EITHER  THE  CUBIC  k-e  OR  THE  LNS 
MODEL 


Accept  and  Exit 


The  figure  shown  above  lists  the  equations  chosen  to  simulate  the  flow.  The  two- 
equation  k-epsilon  model  (RANS)  was  used.  The  earner  phase  equations  were  solved 
using  preconditioning  to  minimize  numerical  diffusion.  The  settings  used  in 
preconditioning  the  equation  are  shown  in  the  following  figure. 


8 


Metacomp  Technologies,  Inc 


FA9550-08-C-0014 


STTR  Phase  II  -Final  Report 
Combustion  Stability  Innovations  for  Liquid  Rocket 


FBemano  Solvers: 

Hrtp  j 

MMmum  dissipation  on: 

♦  LHS  only  v  LHS  and  RHS 

Activate  pressure  acceleration: 

v  y«*  ♦  no  Acceleration  factor. 

{2.0000008+01 

Precond.  Lower  bound  based  on: 

v  Fraction  of  max.  velocity 

♦  Specific  given  velocity 

Transient-Case  precon.  limiting: 

v  Acoustic  CFL-based  Smiting 

♦  No  limiting 

Precond.  turned  off  inibaly: 

v  yes  (give  steps  below) 

♦  no 

Preconditioning  Uscous  Limit  (Set  to  zero  for  Stokes  type  flow) 

0.050 

Precond.  lower  bound  velocity: 

|5.00Q0QGe+01  m/s 

CBck  Here  For  Help  Wth  Preconditioning 

Use  the  LHS  and  RHS"  minimum  dissipation  controls  for  problems  in  which  severe  transients  exist  (e.g  strong  shocks  and 
expansions).  Make  sure  that  for  your  final  solution,  you  have  ramped  DOWN  the  minimum  dissipation  to  0.0.  Use  the  *He^» 
set  CDurant  *  *  button  located  In  the  time  integration  window  to  get  recommended  values  for  *LHS  and  RHS  parameters  and 
to  set  them. 


Accept  and  Exit  | 


A  simple  one-step  chemistry  as  shown  below  in  the  reactions  window  was  employed. 


Reaction  number 


r 


Reaction  forvard  rates  are  calculated  as  follows: 

Rate  -  A*pov(T,N)*pov(PJK)  *exp(-B/  (RoT) ) 

A  ■  Freq.  Fac.  N  Temper.  Exp.  B  -  Act.  Ener. 

P  *  Pressure/ 101 32L . 000000  K  Pressure  Exp.  Unless  you  are  modifying  the 
reactions,  the  coefficient  and  exponent  entries  in  the  reactions  belov  MUST  be 
the  sane.  In  thermal  non- equilibrium  mode,  the  temperature,  T  is  defined  as: 


T  -  pov(Tt,a)  *pov(Tv,b) 

Tt  -  Translational  Temperature 
Tv  -  Tib rational  Temperature 
a  -  Trans,  temp,  expon. 
b  *  Tib.  temp,  expon. 


This  Is  reaction  number  1 

Check  reaction  mass  balance 

Reactants 

Products 

±ld 

+ld 

Coef.  Expon. 

->  -1 

Coef.  Expon. 

[Fo  [Fo  Jh2 

±J  + 

|2  0  jz  0  |H20 

[To  (To  (02 

JJ 

Frequency  factor: 

|50e*13 

sA-1  (mA3/kmoi)^step  reactants -1  (always  SI) 

Temperature  exponent: 

(oo 

Pressure  exponent: 

fo  o— ' 

Activation  energy: 

|0.3e+08 

JTkmoi  or  N  m/kmol  (always  SI) 

Third  Body  Efficiencies  for  Reaction  1 

Accept  and  Exit  | 

WMte  to  FHe 

Help  j  Oose 
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The  options  chosen  to  model  the  dispersed  phase  are  shown  in  the  following  figures. 

p  Dispersed  Species  Properties  x 


Hease  choose  the  dispersed  phase  species  for  which  you  wish  to  enter  the 
properties.  Vou  can  do  this  by  dou  We -clicking  the  left  mouse  button.  If  you 
wish  to  edit  a*  the  species,  then  eftek  on  the  Do  AR  button. 

Note  that  the  drag  force  and  energy  transfer  terms  (when  ap  pic  able)  are 
always  active  In  the  EOF  mode. 

For  unsteady  flows,  when  solving  the  momenUims  and  energy  In  conservation 
form,  please  set  the  Initial  steps  in  non-cons,  form  to  zero 


0 


Solve  momentum  and  energy  m  conservation  form: 
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The  following  figure  lists  the  boundary  conditions  used  in  the  simulation. 
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The  frozen  flow  field  was  initially  setup.  A  heat  source  was  then  applied  to  the  near- 
vicinity  of  the  injectors  to  establish  the  flame.  The  heater  was  removed  after  the  flame 
anchored  itself.  The  time  integration  mode  was  switched  to  unsteady  as  shown  in  the 
following  window  and  calculations  were  continued  until  residuals  reached  a  plateau. 
The  overall  flux  balance  at  this  stage  was  within  1%  of  the  inflow. 
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The  next  set  of  figures  show  the  flow  field  at  different  streamwisc  locations.  The  labels 
on  the  contour  slices  are  x  locations  in  meters.  The  streamwise  velocity  distribution 
below  shows  the  development  of  the  shear  layer  between  the  fuel  and  oxidizer  streams. 


Gaseous  02  at  the  corresponding  planes  can  be  seen  in  the  next  figure.  Almost  50  %  of 
the  O2  fluxes  out  of  the  domain  due  to  lack  of  gaseous  fuel.  The  last  plane  in  the 
following  figure  shows  the  amount  of  O2  leaving  the  test  section.  Lack  of  complete 
combustion  may  be  attributed  to  mixing  of  the  fuel  stream  with  the  inert  coolant  stream. 
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Hydrogen  mass  fractions  (contour  levels  between  0  and  0.1)  at  the  corresponding 
locations  are  shown  below.  About  50  %  of  H2  spills  in  to  the  coolant  stream  and  fluxes 
out  of  the  domain  unreacted. 
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The  next  figure  shows  temperature  contours.  Flame  location  on  the  symmetry  plane  is 
also  visible. 


Temperature  distribution  on  the  wall  is  seen  in  the  next  figure.  The  wall  temperature  of 
almost  700  K  is  reached  towards  the  exit  of  the  domain.  Adiabatic  wall  boundary 
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condition  was  used  at  the  wall  for  this  simulation  in  order  to  estimate  the  maximum 
temperature  at  the  wall.  The  coolant  inflow  boundary  was  subdivided  into  two  sections, 
one  closer  to  the  walls  and  the  other  surrounding  the  three  jets.  The  next  set  of 
simulations  will  alter  the  inflow  velocity  in  the  region  next  to  the  wall  to  estimate  the 
effect  of  coolant  flow  rate  on  maximum  temperature  at  the  w  all. 

A  second  case  at  higher  chamber  pressure  produced  higher  wall  temperatures  as  shown 
below.  The  maximum  temperature  increased  to  1000  K  at  the  wall.  In  this  case  flow 
remained  unsteady  and  only  an 


1000  K 

t 


instantaneous  snapshot  is  shown. 

Conclusion 

The  capability  in  CFD++  to  resolve  acoustic  oscillations  and  their  coupling  with  reaction 
zones  was  improved.  CFD++  can  now  resolve  acoustic  oscillations  in  combustion 
chambers.  Acoustic  modes  were  obtained  for  AFRL  combustor  stability  test  rig  under 
quiescent  conditions  to  study  the  impact  of  wave  guide  on  the  acoustic  response  of  the 
fucl/oxidizer  jets.  Steady  state  calculations  of  the  test  facility  are  continuing  to  proceed 
in  order  to  asses  heat  load  on  the  wave  guide  for  hot  fire  eases. 
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Effort  and  Accomplishments  at  Penn  State  University  and  Georgia  Tech 

Cry  ogenic  Fluid  Dynamic  Response  of  Swirl  Injector  to  External  Forcing  at 
Supercritical  Conditions 

This  task  extends  our  previous  study  on  the  eryogenie  fluid  dynamics  of  swirl  injeetors  at 
supercritical  conditions  to  investigate  the  dynamic  response  of  the  injector  to  external 
excitations.  The  forcing  is  imposed  through  periodic  oscillations  of  the  mass  flow  rate  at 
the  injeetor  tangential  inlet  over  a  wide  range  of  frequencies  relevant  to  the  intrinsic  flow 
instabilities  of  the  injector. 

The  dynamic  response  of  swirling  jets  to  external  excitations  has  been  investigated  by 
Panda  and  McLaughlin1  using  flow  visualization  techniques.  When  the  jet  was  excited  by 
acoustie  oscillations  with  discrete  frequencies  in  the  upstream  region,  organized 
structures  evolved  from  the  originally  weak,  irregular,  and  large-seale  structures  under 
conditions  without  forcing.  Cereeedo  et  air  found  that  when  a  eo-flow  jet  was  forecd  at 
its  natural  or  sub-harmonie  frequencies,  the  turbulence  properties  were  significantly 
ehanged.  However,  at  the  foreing  frequency  twiee  of  the  fundamental  value,  only  modest 
influence  was  observed.  In  the  experimental  study  of  Gallairc  et  al},  the  vortex 
breakdown  phenomenon  was  found  unaffected  by  azimuthal  foreing  at  the  nozzle 
periphery  since  this  proeess  was  governed  by  the  dynamics  in  the  core  region.  Bazarov 
and  Yang4  developed  a  generalized  linear  theory  on  injector  dynamics.  Richardson  et  al / 
studied  the  dynamic  response  of  a  simple  swirl  injector  using  a  boundary  clement  method 
and  compared  the  results  with  linear  theory. 

So  far,  very  limited  effort  has  been  applied  to  study  eryogenie  fluid  dynamics  in  a  swirl 
injeetor  using  high-fidelity  modeling  and  simulation  techniques,  especially  at  conditions 
representative  of  liquid  rocket  operation.  The  present  task  is  to  develop  a  LES  (large 
eddy  simulation)  based  numerical  framework  within  which  various  physical  processes 
and  flow  parameters  dictating  the  dynamic  response  of  a  cryogenic-fluid  swirl  injeetor 
can  be  studied  systematically.  The  theoretical  formulation  and  numerical  approach 
follow  those  reported  in  Zong  and  Y ang*  The  injeetor  geometry  and  boundary  conditions 
are  identical  to  the  baseline  ease  of  our  previous  study.*  The  mean  mass  flow  rate  is 
0.15kg/s.  External  forcing  is  imposed  by  pulsating  the  mass  flow  rate  at  selected 
frequencies  at  the  tangential  inlet  of  the  injector: 

(1) 

where  ^oand  If  denote  the  mean  mass  flow  rate  and  foreing  frequency,  respectively, 
and  a  the  oscillation  magnitude,  which  is  fixed  to  be  15%.  The  foreing  frequency  covers 
a  range  of  0.55-14  kHz,  in  order  to  eapture  the  broadband  nature  of  the  injeetor  flow 
dynamics.  Turbulenee  with  an  intensity  of  8%  of  the  mean  flow  quantity  is  provided  by 
super-imposing  broadband  noise  onto  the  instantaneous  mass  flow  rate. 

Figures  1-2  show  the  temporal  evolution  of  the  temperature  fields  at  foreing  frequencies 
of  14  kHz  and  3.2  kHz,  respectively.  The  14  kHz  foreing  does  not  show  any  significant 
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influence  on  the  flow  dynamics.  The  wave  structure  resembles  that  of  the  unforced  ease. 
Small  disturbances  initiated  in  the  upstream  region  grow  as  they  are  conveetcd 
downstream.  As  noted  in  our  previous  study*,  the  wave  structure  features  three- 
dimensional  hydrodynamie  instability  and  propagates  in  a  form  similar  to  the  shallow- 
water  wave  for  an  incompressible  flow.  The  high  frequency  fluctuation  decays  rapidly 
and  disappears  as  the  LOX  film  moves  downstream  due  to  viscous  dissipation. 
Consequently,  foreing  at  this  frequency  does  not  affect  the  flow  field  significantly.  For 
the  3.2  kHz  foreing  ease,  the  wave  structures  on  the  LOX  film  lose  their  identities  shortly 
after  their  appearance.  They  breakup  and  merge  into  larger  billows.  The  frequency 
spectra  of  pressure  fluctuation  provide  more  quantitative  insight  into  the  various  types  of 
instability  waves  in  the  injector  flow.  Figures  3-5  indicate  that  the  fluctuation  magnitude 
at  the  foreing  frequency  stands  up  in  the  bulk  of  the  LOX  film.  In  addition,  forcing 
makes  the  power  of  other  frequency  components  smaller  and  the  frequency  content  is 
reduced  greatly,  a  situation  commonly  known  as  the  frequency-locking  phenomenon. 
Even  the  originally  dominant  fluctuation  may  lose  most  of  its  energy  at  some  foreing 
frequencies. 

Figures  6-9  show  the  variations  of  the  mass  flow  rate  at  the  injector  inlet  and  exit  for  the 
eases  with  3.2  kHz  foreing  and  without  oscillation.  Foreing  causes  the  mass  flow  rate 
and  film  thiekness  to  vary  in  a  more  regular  manner  and  the  fluctuation  magnitude  of  the 
spreading  angle  becomes  larger.  Figures  10-11  show  the  time-averaged  film  thickness 
and  spreading  angle  as  a  function  of  foreing  frequency.  The  mean  film  thickness  and 
spreading  angle  do  not  change  much,  due  to  conserv  ation  of  mass  and  momentum. 

The  injection  dynamics  ean  be  quantified  globally  in  terms  of  the  mass  transfer 
function,  defined  as: 

I K»=§  (2) 

where  ma  is  the  Fourier  component  of  the  mass  flow-rate  fluctuation  at  the  forcing 
frequency.  The  subscripts  ex  and  in  denote  the  injeetor  entrance  and  exit,  respectively. 
Figure  12  shows  the  magnitude  of  the  injeetor  mass  transfer  function  as  a  function  of  the 
foreing  frequency.  This  quantity  reaches  its  maximum  around  1.9  kHz,  corresponding  to 
the  dominant  frequency  of  the  intrinsic  injeetor  flow  instability  under  the  condition 
without  external  foreing. 
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Figure  1 .  Temporal  evolution  of  temperature  field,  external  forcing  at  14  kHz. 
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Figure  2.  Temporal  evolution  of  temperature  field,  external  forcing  at  3.2  kHz 
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Figure  3.  Power  Spectral  densities  of  pressure  fluctuations  at  different  locations  inside 

injector,  external  forcing  at  1.9  kHz 
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Figure  4.  Power  Spectral  densities  of  pressure  fluctuations  at  different  locations  inside 

injector,  external  forcing  at  2.1  kHz 
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Figure  5.  Power  Speetral  densities  of  pressure  fluctuations  at  different  locations  inside 

injector,  external  foreing  at  3.2  kHz 
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Figure  6.  Mass  flow  rates  at  injeetor  inlet  and  outlet  without  external  foreing 
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Figure  7.  Mass  flow  rates  at  injeetor  inlet  and  outlet  with  external  foreing  at  3.2 

kHz 
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Figure  8.  Comparison  of  spreading  angle;  free  oscillation  and  3.2  kHz  forcing. 
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Figure  9.  Comparison  of  injector  exit  film  thickness;  free  oscillation  and  3.2  kHz 

forcing. 


Figure  10.  Variation  of  spreading  angle  with  forcing  frequency 
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Figure  1 1 .  Variation  of  film  thickness  at  injector  exit  with  forcing  frequency 


Figure  1 2.  Frequency  spectrum  of  magnitude  of  injector  mass  transfer  function 
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